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Families of explicit solutions are found to a nonlinear Black-Scholes 
equation which incorporates the feedback-effect of a large trader in case 
of market illiquidity. The typical solution of these families will have a 
payoff which approximates a strangle. These solutions were used to test 
numerical schemes for solving a nonlinear Black-Scholes equation. 



< 

c^ ' Key words and phrases: Black - Scholes model, illiquidity, nonlinearity. 

explicit solutions 
AMS classification: 35K55, 22E60, 34A05 

> 

o 
o 



X 



^e-mail: Ljudniila.Bordag@ide.hh.se 
^e-mail: chniakova@math.tu-cottbus.de 



1 Introduction 

Standard option pricing theory uses a number of basic assumptions including the 
assumptions of symmetric information, of complete and frictionless markets, as 
well as the assumption that all participants act as price takers. Recently a series 
of papers appeared in which one or more of these assumptions have been relaxed; 
[T^ . [T7j . [T], [IH], IS] and j2] are representative examples of this work. The 
turbulence on financial markets such as the events surrounding the collapse of 
LTCM in 1998 have made market liquidity an issue of high concern for investors 
and risk managers and have triggered a lot of academic research; see for instance 
[T^ . [^, P3]. In illiquid markets an attempt to buy/sell a large amount of an 
asset will affect its price so that the assumption that investors act act as price 
takers cannot be maintained. 

The purpose of this paper is to investigate the evaluation of an option hedge- 
cost under relaxation of the price-taking assumption. For our analysis we use 
the framework proposed by Frey in j3], |3]. He developed a model of market 
illiquidity describing the asset price dynamics which result if a large trader chooses 
a given stock-trading strategy {at)t- The resulting stock-price dynamics have the 
following natural property: if the large trader buys (sells) stock, i.e., if daj > 
(dtti < 0) the stock price rises (falls). If the position of the large trader is 
unchanged, the stock price St follows standard geometric Brownian motion with 
constant volatility a. Formally, Frey models stock price dynamics by the following 
stochastic differential equation 

dSt = aSt^dWt + pSt^dat, (1.1) 

where Wt is a standard Brownian motion and St- denotes the left limit \im.s^t,s<t St. 
In ()1.1|1 p is the market illiquidity parameter with < p. The value l/{pSt) is 
called depth of the market at time t. Note that in the model (jl.lj) the parame- 
ter p is a characteristic of the market and does not depend on the payoff of the 
hedged derivatives. If p — > then ()1.1|) reduces to the Black-Scholes model. We 
concentrate our investigations on the nontrivial case p 7^ 0. 
Consider the problem of hedging a terminal-value claim with maturity T and 
payoff h{S) in the model (jl.lj) . As shown in |3], |lj, the feedback-effect leads to 
a nonlinear version of the Black-Scholes partial differential equation for a hedge 
cost u{S,t) of the claim, 

with terminal condition u{S,T) = h{S). The variable S denotes the price of the 
underlying asset and t is the time variable. The equation above is studied for the 
variables S and T in the intervals 

S>0, te [0,T], T>0. (1.3) 



Similar equations in related models were obtained by a number of authors see for 
instance P, 0, HSj, [HI, [H, [HI, [H]. 

Frey and co-authors , |lj , [Hj studied equation ()1.2|) under constraint and did some 
numerical simulations. Our goal is to investigate this equation using analytical 
methods. 

We study the model equation p.2p using methods of Lie group theory in Section|21 
Using the symmetry group we reduce the partial differential equation p.2j) to an 
ordinary differential equation in Sectional We obtain nontrivial explicit solutions 
for this case. We prove that the explicit solutions approximates strangles with 
corresponding payoffs (see Section IH). Further, in Section El we study different 
properties of the obtained solutions. The existence of nontrivial explicit solutions 
allows us to test different numerical methods usually used to calculate hedge-costs 
of derivatives. The best results are achieved by the completely implicit method. 
The validated numerical scheme was used to calculate option hedge-costs in case 
of calls and bull-price-spreads. 

2 Lie group symmetries 

In this section we study the symmetry properties of equation ()1.2j) and obtain the 

complete description of the corresponding Lie algebra, the associated Lie group 

and a list of functionally independent invariants. 

Let us study the nonlinear part of this equation. The denominator in the second 

term of this equation will be equal to zero if the function u{S,t) satisfies the 

equation 

1 - pSuss = 0. (2.4) 

The solution of this equation is a function Uo{S, t), 

Uo{S,t) = -S\nS + Sciit) + C2{t), p^O, (2.5) 

P 

where the functions ci(t) and C2(t) are arbitrary functions of the variable t. From 
now on we assume that the denominator in the second term of equation ()1.2|1 is 
not identically zero, i.e., the function u{S,t) is not equal to the function Uo{S,t) 
()2.5|) except in a discrete set of points. 

We introduce the necessary notations connected with the Lie group theory. Be- 
sides the classical work (TU] our notations follow |H] and, especially with respect 
to the invariants, to Ovsiannikov ^21 and Olver TT]. We introduce the two- 
dimensional space X of independent variables {S, t) G X and a one-dimensional 
space of the dependent variables u ^ U. Then we consider the space f/(i) of the 
first derivatives of the variable m on S" and t, i.e., {us, Ut) G t/(i). Analogously we 
introduce the space U(2) of the second order derivatives {uss,ust,Utt) G t/(2). Let 
M = X xU he the Cartesian product of pairs (x, u) with x = (5, t) G X, -u G U. 



The second order jet bundle M^"^^ of the base space M has the form 

M(2) = X X f/ X f/(i) X t/(2). (2.6) 

We label the coordinates in the space M^^^ hy w = {S,t,u,us,Ut,uss,ust,Utt) G 
M^^\ The second order jet bundle M^"^^ has a natural contact structure (see |12j . 
[TT] . [in], jH], 01 )• Our differential equation ()1.2p is of order two and in the 
context of the second order jet bundle M^^^ it should be seen as an algebraic 
equation in M^'^\ We introduce the following notation, 

A{S, t, u, Us, ut, uss, ust, utt) =ut-\ -—-- r^. (2.7) 

Equation p.2j] is then equivalent to the relation 

A{w)=0, weM^^l (2.8) 

We identify this algebraic equation with its solution manifold La defined by 

La = {u;G M(2)|A(w) =0} cM(2). (2.9) 

We consider an action of Lie-point groups on our differential equation and its 
solutions. We are interested in the group Diff(M(^'') compatible with the contact 
structure of M^'^\ We denote the corresponding algebra by Viff^M^"^^). The 
symmetry group Ga of A is defined by 

GA = {^eDiff(M(2))| g: La^La}. (2.10) 

Theorem 2.1 The differential equation Hl.^) admits a nontrivial four dimen- 
sional Lie algebra spanned by generators 

d d d d d 
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Proof. Let us consider a Lie-point vector field on M, whose elements are repre- 
sented by 

C\ C\ C\ 

F = e(^,t,M)— + r(5,t,M)- + 0(5,t,u)— , (2.11) 

where ^{S, t, u), t{S, t, u) and 0(5', t, u) are smooth functions of their arguments, 
V G Viff{M). Assume there exists an infinitesimal generator of an action 
g e Ga- The infinitesimal generators of these transformations form an algebra 
Dif f ^{M). A Lie group of transformations acting on the base space M induces 
transformations on M*^^^ The corresponding algebra Viff^^M^"^^) will be com- 
posed of the vector fields 

pr^'^V = ^{S,t,u)-^ + T{S,t,u)^^+<P{S,t,u)-^ 

+ ^S(S,t,u)j— + (P\S,t,u)j— (2.12) 

dus dut 

+ 0^5(5', t,M)-^ + 0^*(5,t,n)-^ + 0«(5,t,n)-^ 
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where pr^'^'^V is the second prolongation of the vector filed V. Here the smooth 
functions 4>^{S, t, u), 0*(S', t, u), (p^^{S, t, u), (p^^{S, t, u) and 0**(S', t, u) are uniquely 
defined by the functions ^(5*, t, u), t{S, t, u) and 0(5', t, u) using the prolongation 
procedure (see [12], [H!, [lEl, 0, 0)- 

For our calculations we will use the explicit form of the coefficients 0*(S', t, u) and 
(f>^^{S, t, u) only because of the special structure of equation ()1.2|) . The coefficient 
0*(5', t, u) can be defined by the formula 

0*(5', t, m) = 04 + Mt0„ - U5^i - UsUtiu - Utn - {UtfTu (2.13) 

and the coefficient (p^^{S, t, u) by the expression 

0^^(5,t,M) = 055 + 2m505« + W550« (2.14) 

+ {Usf(t)uu - '^Ussis - UsisS - 2(m5)^^5„ 

- 3M5M55^„ - {Usfiuu - 2ustrs - UtTss 

- 2usUtTsu - {UtUss + 2usUst)Tu - (usYutTuu, 

where the subscripts of C,, t, denotes corresponding partial derivatives. The 
symmetry algebra Viff^^M^"^^) of the second order differential equation A = 
can be found as a solution of the determining equation 

pr(2)l/(A) = {mod{A = 0)), (2.15) 

i.e., the equation ()2.15j) should be satisfied on the solution manifold La- It is 
easy to prove that equation ()2.15j) has the following solutions, 

ou ou 

where Vi G T>iff^{M),i = 1,2,3,4. The commutative relations are 

[Vi, V2] = [V,, Vs] = [V2, V,] = [V2, V,] = [Vs, V,] = 0, 

[Vi,V4 = -V^. (2.17) 

The vector fields Vi, i = 1, 2, 3, 4 span a four dimensional solvable Lie algebra. 
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An element of the algebra T>iff^{M) can be represented as a linear combination 
of the vector fields given by formulas ()2.16|) 

d d d 

V = aiVi + 02^ + OsVs + 04^4 = ^a{S,t,u)— + Ta{S,t,u)— + 4>a{S,t,u) — , 



where 

^a{S, t, U) = aiS, Ta{S, t, U) = 02, (pa{S, t, u) = ttiM + a^S + a^ 

with arbitrary constants Oi, 02, 03, 04. 

Every element V of the algebra Viff^{M) is an infinitesimal generator of an 

action g G Ga- Using the Lie equations we prove the following theorem. 

Theorem 2.2 The action of the symmetry group G/\ of U.^) is given by I12.20\} ~ 

Proof. To find the transformations of the Lie group Ga associated with the 
generators ()2.16|1 we just integrate the system of ordinary differential equations, 
the so-called Lie equations, 

— =^a{S,t,U), — = Ta{S,t,U), — = (j)a{S ,t, u), (2.18) 

with initial conditions 

S\,=o = S, t|e=o = t, m|,=o = u, (2.19) 

where e is the group parameter. Here the variables S, t and u denote the values 
S, t, u after a symmetry transformation. The solutions to the system of ordinary 
differential equations ()2.18|) with initial conditions ()2.19p have the form 

S = Se^'', eG (-00,00), (2.20) 

i = t + a2e, (2.21) 

u = ue"i^ + a35ee"^^ + — (e"i^-l), ai^O (2.22) 

ai 

u = u + a^Se + a^e, ai = 0. (2.23) 

The equations ()2.20|) - ()2.23|) represent the action of the four parametric symmetry 
group Ga- 
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We will use this symmetry group to construct invariant solutions to equation 

(ll.2|l . In detail the method of construction of invariant solutions is given in the 

book |12j and in the third chapter of the book [llj. A lot of examples are given 

in the books [HI, 0, 0- 

To obtain the invariants of the symmetry group Ga we can use a shortcut because 

of the very simple structure of the Lie algebra found. 

We exclude e from the equations ()2.2()|1 - ()2.23|1 . Two functionally independent 

invariants can be taken in the form 

invi = ait — 02 In S*, (2.24) 

inv2 = ai— — a^lnS + —,S > 0. (2.25) 



The functions (I2.24|) - ()2.25|) are not defined at the point S = and, although the 
model equation ()1.2|) is defined at that point, we will exclude S" = in all further 
investigations. 

We remark that the form of these invariants is not unique. Each function of 
invariants ()2.24j) . ()2.25|1 will be an invariant. Especially we can multiply each of 
the invariants by a constant because any constant is a trivial invariant of the group 
Ga- But it is possible to obtain just two nontrivial functionally independent 
invariants which we take in the form ()2.24|) . ()2.25p . The invariants can be used 
as new independent and dependent variables. 

3 Scaling variables 

Using the symmetry group Ga found in the preceding section we reduce equation 
()1.2|1 to an ordinary differential equation and define families of invariant solutions. 

Theorem 3.1 Up to the group transformations given by i2.2(J\} - r2.2^} all non- 
trivial Lie invariant solutions to equation ^l.S^) depend on the scaling variables 
z,v{z) and the relations 

z = \nS-6t, 6^0, (3.26) 

u{S,t) = -Sv{z), (3.27) 

where 6 is an arbitrary constant, hold. 

Proof. We can reduce the partial differential equation p.2|) to an ordinary 
differential equation for the function v{z) if we change the variables u,S,t for 
z = (f){S,t,u) and v = tp{S,t,u). This substitution leads to invariant solutions 
to equation ()1.2|) if (f){S,t,u) and ilj{S,t,u) are some invariants of the symmetry 
group Ga- In the previous section we found just two invariants, hence all invariant 
solutions except for trivial ones will arise after the substitutions ()3.26|) - ()3.27|) . 
The trivial solutions we can obtain if we assume u = const., u = u{t) and 
u = u{S). 

We remark that we take as a new independent variable the first invariant ()2.24j) 
of the symmetry group Ga and as the dependent variable the nontrivial part 
of the second invariant, this allows us to simplify the calculations. In this way 
we do not lose any solutions because the found invariant solutions can be later 
transformed by the rule of thumb given by ()2.22|1 - ()2.23|1 . 
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The equation for the function v{z) has the form 

2 

v,{l + p {v, + v,,)f -^Avz + V,,) = 0. (3.28) 

zo 



The Lie group of symmetries for this equation can be found in the same way as 
described in previous Section for equation ()1.2|) . 

Theorem 3.2 UUf The equation \3.2^) admits a two dimensional Abelian Lie 
algebra spanned by two generators 

Proof. This theorem was proved in a more general case by S.Lie in (llOj. Also 
it can be verified by a straightforward calculation. 

n 

Equation ()3.28|) allows a two-dimensional Lie group associated with the Lie al- 
gebra spanned by the generators ()3.29|) . As a consequence equation ()3.28p is 
completely integrable. Hence the most general form of the solution of ()3.28|) is a 
two parametric family of congruent curves. To obtain a two parametric family of 
solutions to equation ()3.28p we can subsequently use the two generators ()3.29|) in 
arbitrary order. Both ways will lead to the same family of solutions independent 
on the order. To obtain a solution we must perform two integrations and this 
procedure is not always possible in closed form. However, in view of the theo- 
rem 13.21 we do not have any other possibility to solve equation ()3.28|1 in a more 
convenient way. 

In the next Section we put constraints on the constant 6 in fj3.26|) in order to 
integrate the arising equations in an exact form. Consequently we restricted 
ourselves and do not obtain the most general form for the family of solutions. 



4 Families of invariant solutions 

Theorem 4.1 The equation i3.2^} can be reduced by the substitution v{z) 
y{z) to the set of equations 



y{z) = 0, y{z) = -\-l±\-\, (4.30) 




dy ^ ( f 2 , y ^^ ^ I 1 / ^^ / ^^ 



where 5 is an arbitrary constant. 

The complete set of solutions to equation \3.2^) coincides with the union of so- 
lutions to these equations. 



Proof. First we look for the solutions of the type v{z)z = const. From straight- 
forward calculations we obtain that equations ()4.30|) hold. The corresponding 
solutions to equation ()3.28|) have the form 

viz) = ci, (4.32) 

where ci is an arbitrary constant, and 

1 / [^\ 

v(z) = — 1± a/— 2 + const. (4.33) 

p \ V 2()y 

We assume now that y{z) ^ const., i.e., Vz{z) ^ const, and use the operator U2, 
()3.29|1 . first to introduce the new dependent variable y{z) = Vz{z) in equation 

(IS2H1). 

We obtain a first order differential equation for the function y{z), 

2 , ^y- f 2 , y ^M , /^ 2 , 2^_ , 2(5-ct2 



The equation ()4.34j) is quadratic in the highest derivative and it can be repre- 
sented as a product of two differential equations ()4.3H1 . 

We reduced equation ()4.34j) to a product of two equations ()4.31|1 and in this way 
we could have lost some of the solutions. Let us now study the discriminant 
curve for equation ()4.34|) . We denote by F^y^, y, z) the left hand side of equation 
dlSl, i.e., 

F(y., y. ., ^ yl ^ 2^ {y^ ^ ^ _ -^) ^ (,^ ^ 1 , , ?il^) . ,4,5) 

The discriminant curve is a set of points fulfilling the conditions, 

F{y,,y,z) = Q, (4.36) 

^%^ = 0. (4.37) 

dyz 

Along this curve the conditions of the theorem on an implicit function are not 
satisfied and in these points the obtained solutions may be not unique. It is easy 
to prove that the system of equations ()4.37p - ()4.36|) has a unique solution, 

yexcep{z) = -, (4.38) 

for the special value of the constant 6 

6 = aVS (4.39) 

only. The corresponding solution of equation ()3.28|) has the form 

2; 
Vexcepiz) = - + COUSt (4.40) 

P 

and it coincides with one of the solutions ()4.33|) for S = a'^/8. 
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Theorem 4.2 The explicit invariant solutions to equation M.^) . defined on the 
region S* > 0, t G [0,T], T > are given by \4-4i\ ) A4-4^) '^'^^ \4-46\) - Other 
solutions of this type can be obtained using the transformations of the symmetry 
group Ga represented by \2.2^^) - r2.2!^) . 

Proof. To obtain the invariant solutions we should solve the equations listed in 
the theorem 14. 1[ It is trivial to solve the first two of them. 
The relations ()4.H()j) have the following solutions 

u{S,t) = Sci, (4.41) 

and 



u{S, t) = p-M 1 ± W — (S In S - 5St) + SdQ, (4.42) 

where 5 and do are arbitrary constants. 

To integrate the equations ()4.H4|1 we use the second operator f/i, ()H.29|1 . and 

separate variables 



y cr^ \ 1 / cr^ / (T^ 



-1 



We denote the integral on the left hand side by 



2 . y ^^ ^ I 1 / ^^ / ^^ 



-1 



Straightforward calculations lead to the following form for the function /(y), 



^(^) = 2i3^ -2^^M^^2V^ + ^^^M^^te^-^^ 



0"^ , / _ / ^ , i / Cr 




f ^f ^ 1 



2 / / T ^2 



T V2^"^^'"^M^^^2VV 



where the variables ^ and y are connected by 



a 



2 



- y = r • (4.43) 



Spa 
Now let us chose 6 = cr^/8, i.e., in a way that substitution ()3.2(i|l takes the form 

2 

z = \nS-^t. (4.44) 

8 
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We obtain an explicit representation for the function y{z) which solves equation 
()4.34|) . The solutions are given by 



/ 



yi^) 



V 



24/3g. 

( (m + eiA/m2"+4e3^] j 



1/3 



+ 



( (m + eWm'^ + Ae^^lA j 



1/3 



\ 



24/3( 



(4.45) 

with an arbitrary constant m and Ci = ±1. 

Thereafter we integrate the equation Vz = y{z) and obtain a family of solutions 
to equation ()3.28p . 



/ 



Viz] 



ei\J ni^ + 4e 



_2-4/3g-. j (^ _^ ei\lm? + 4eiM j 



1/3 



In rri + t2\J m? + 4e 



1/3 



^^ + ^2\j "1^ + 4e^ 



1/3N 



+ ci 



where d and m 7^ are arbitrary constants. The case m = corresponds to the 
solution ()4.33|1 . The parameters e\ and e^ take values ei = ±1 , £2 = ±1 and can 
be chosen independently. The solutions do not depend on the value of e\. 
The corresponding family of solutions to equation ()1.2j) will take the form 



0" 



uiS.t) = p-^SlnS St 

8p 



-2-^Vexp(4^^ f {m + eiJm^ + 4S^/^exp(-^ 



4\ 1/3 



4\ 1/3 



.2-4/3-1 



p exp 



ah 



-m 



+ £i\ rn? + 453/2 exp 






-p ^S\n\ m + e2 W m^ + 45*3/2 exp 






1/3 



(4.46) 



-m + 62 W m? + 45*3/2 exp 



3ff2t~ 
~[6 ) 



1/3N 



+ 5(ii + 62, 
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where di,d2 are arbitrary constants. In formula ()4.46|) we assume that the ar- 
bitrary parameter m is non equal to zero. In case m = this solution can be 



reduced to one of the solutions ()4.42|1 . 

Let us compare this solution with the solutions to equation ()4.H4|1 which were 
obtained in the case 5 = a'^/8. The functions y{z) in the family ()4.45|1 are even 
functions of the constant m. For m = we obtain 

y{z) = --. (4.47) 

P 

The solution ()4.47|1 leads to the described solutions ()4.33|1 and ()4.42|1 with upper 
sign and with 6 = a'^/2. 
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5 Properties of invariant solutions 

The solutions ()4.4(jj) depend on the parameter p in a very simple way: all solutions 
of this family have the factor 1/p in front of the whole expression. This parameter, 
which is the measure of the influence of the large trader on the market, is a 
constant < p and cannot be equal to zero for the large trader. This means that 
each solution of this family does completely blow up at p — > 0. Consequently 
these solutions have no linear analogies. 

If we denote by u{S,t) = pu{S,t)p we obtain for the function u{S,t) following 
equation 

This means that the solutions ()4.4(jj) multiplied by p are solutions to equation 
()5.48|) . If we obtain any solutions to ()5.48|) for any fixed boundary conditions, we 
obtain the corresponding solutions to ()4.46|) with boundary conditions divided 
by p if we divide the found solutions by p as well. In other words, the solutions 
()4.4fi|l strongly reflect to the nonlinearity in equation ()1.2j) . 
Let us study the analytical properties of the solution ()4.4(ij) and the corresponding 
payoff. In Figure^ we represent graphically the solution u{S,t) ()4.46p for small 
values of the variables S, t. 

Let us represent the payoff of a strangle as a sum of Kp European puts with an 
exercise price Ep and Kc European calls with an exercise price Ec which have 
the same expiry date T. We can choose in an appropriate way the parameters 
m,di,d2 of the explicit solution ()4.46|) such that this solution approximates the 
payoff of a strangle UstrangUS, T), 

UstrangieiS, T) = Kp max {Ep -S,0) + Kc max {S - Ec, 0), Ep < Ec- (5.49) 

12 



This is shown in Figure |21 

Let us now investigate the asymptotic properties of solutions (I4.46|) for S* ^ 
and for S ^ oo. The asymptotic behaviour of the function ()4.46|1 for S ^ can 
be described as follows, 

pu{S,t) ~ (m^)^/^exp(^tj +51n5 (5.50) 

The main term in formula ()5.50|) depends on the time and on the constant m. 
We can choose m to model payoff properties. From this decomposition it follows 
immediately that for all solutions from this family the denominator in equation 
(ll.2|l vanishes in the point 5 = 0. In order to avoid this singularity we exclude 
the point S = from the intervals where the numerical investigations are done. 

The main term of the asymptotic expansion of u{S, t) for S -^ oo, 



exp('^')m2 + 0(5-^/^), S ^ oo, (5.51) 



is equal to 35* In 5*. This term is independent of any integration constant or time. 
Hence all solutions in this family have the same asymptotic behaviour for 5* — > cxd. 
From the financial point of view it is important to study the dependencies of the 
obtained solutions on different parameters, for instance, on time, on volatility or 
on the price of the underlying asset, etc. In this way we get information about 
the sensitivity of our product with respect to a change of one of these parameters. 
Using the explicit formula for the solutions ()4.4fj|l it is easy to represent these 
dependencies graphically, see Figures El -El The time dependence of the solutions 
()4.46|) is very weak but still present as we can see on Figure JHl 
The obtained family of solutions ()4.46|) can be used as a benchmark for testing 
of numerical methods. We suggest the following procedure. We use the solutions 
()4.4fj|l with boundary conditions which we can obtain just by fixing the time to 
test numerical methods. These boundary conditions are smooth. Then we take 
one of the numerical methods and try to reproduce the analytical solution. In 
this manner we can check on each time step the reached accuracy and adjust 
the parameters of the grid and the numerical scheme. Thereafter we can be sure 
that for all smooth boundary conditions of the same type as studied we obtain 
numerical solutions with nearly the same accuracy. 

Now if we apply this method to boundary conditions with worse properties we can 
be sure that it works at least in the case of an approximation of these boundary 
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conditions by very close but smooth ones. In case of an European call option we 
have a continuous payoff function Ui{S, T), 

ui{S,T)=max{S-E,0), (5.52) 

where E is the exercise price and T is the expiry date, which is not differentiable 
in S* = ii^. We can make it smooth by just replacing the payoff Ui{S,T) in the 
neighbourhood of the exercise price E. Usually one takes as such smooth func- 
tion a solution of the linear Black-Scholes formula ()5.5H|1 for t ~ T. Then we can 
compare the results of numerical calculations in both cases. If they do not have 
any significant difference we can use the same method also in case of continuous 
boundary conditions and relax the condition of smoothness. 

As a first example we take an explicit method for a numerical solution of equation 
(ll.2|l . This method can be used to find numerically solutions to the linear Black- 
Scholes model 

2 

Ut + —S'^uss + rSus - ru = 0, (5.53) 

where r is the interest rate. It gives proper results for the special relation be- 
tween AS*^ and At, where by AS, At we denote correspondingly the mesh sizes of 
the discretization of 5* and t intervals. We applied this method to the nonlinear 
equation p.2j) . We proved that in all studied cases the explicit method diverges 
independently from the chosen relation between AS^ and At. It follows that the 
explicit method is not reasonable in this nonlinear case. 

Another way to solve equation p.2j) numerically is to use the completely implicit 
method. For the linear Black-Scholes model ()5.53j) it gives proper results for 
arbitrary relations between AS"^ and At. For a nonlinear equation this method 
leads to a system of nonlinear algebraic or transcendental equations. An attempt 
to solve such a system can easily exceed the possibilities of a modern computer 
due to the very fast with growing grid size. We used this method to reproduce the 
explicit solutions ()4.46|) with appropriate accuracy. Thereby we used equidistant 
grids with 16, 28 and 42 space nodes and with 15 and 30 time levels. We reached 
the relative accuracy of order of 0.2% . 

Then we used this completely implicit method to calculate the value of derivatives 
governed by equation ()1.2p with usual payoff functions. 

Let us describe shortly the system of difference equations which we used. It was 
obtained by replacing the derivatives in the t and S directions in the following 
way. 



du _ u{Si,tj+i) -u{Si,tj) 
dt 



T 
^2 



+ 0(r), 

+ 0(/i'), (5.54) 



d^u _ u{Si+i,tj) - 2u{Si, tj) + ujSi^i, tj) ^ ^^^2 
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where r = At is the time step and h = AS the space step. For each fixed j we 
obtain a system of Ns — 1 equations 

— - p[Ui^ij - 2uij + Ui+ij) T "c^ ~ P{Ui^ij - ^Uij + Ui+ij^ 



4 V^,- "^^ '-'' '^ ' '^^^V 4 VS 



H — {ui-ij - "^Uij + Ui+ij) = 0, i = 2, Ns, j = Nt, 1, (5.55) 

o 

for the internal points, where Nt is the number of time layers and Ns + 1 is the 
number of grid nodes in space direction. In this case we used the final conditions, 
i.e., the knowledge of the values u{S,T) and calculated the values for u{S,t) 
backwards to t = 0. In the system ()5.55|) the values on the layer j + 1 are known 
and the values on the layer j are unknown functions. On the boundaries 5*1 and 
Sns+1 the values Uij are defined for each fixed j by the function UboundiS, t) in 
accordance with the used boundary conditions. The complete system of difference 
equations has the form 



TT - P{Ubound{Sl, tj) - 2U2J + U3j) 



4 V5. 



{Ubound{Sl,tj) - 2U2j + U3j) = 0, (5.56) 



,2 



Uij+l — Uij f h" n , 

— - p{ui^ij - 2uij + Ui+ij^ 



ra'h 



2U2 



-\ — {ui-ij - 2uij + Ui+ij) = 0, i = 3,Ns-l, (5.57) 



'^TV.i+l — UNsJ I h 



2 

— p{uNs-lj — '^UNsJ + Ubound{SNs+li tj)) 



H {uNs-lj - '^UNsj + UboundiS Ns+lj,tj)) = (5.58) 



with j = Nt, 1. 

In the works [3] and |3] it was proved that the hedge-cost of the claim u{S,t) 
increases monotonously with growing p, i.e., with growing infiuence of a large 
trader. We prove this dependence numerically. We take as boundary conditions 
UboundiS, T) = ui{S,T) ()5.52|1 . i.e., the boundary conditions which correspond to 
one European call option. We calculate the values u{S, t = 0) for various values of 
p. In Figures EI-IHl we can see that with the growing value p the option hedge-cost 
also grows monotonically. It completely corresponds to the functional behaviour 
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obtained in the work j3]. 

Now we compare the option hedge-cost predicted by the hnear Black-Scholes 
model ()5.53|) and by the nonhnear model ()1.2p . 

At first we find numerically the value of the hedge-cost for the derivative u^lS, 0) 
defined by equation p.2|) with the payoff ^3(6', T) which corresponds to three 
European call options. The payoff function for K European call options is given 

by 

UKiS,T) = Kmax{S-E,0), (5.59) 

where we will use K = 3, 5, 8. Then we find numerically the value of the hedge- 
cost for the derivative u^^S, 0) defined by equation p.2j) with the payoff u^i^S, T) 
()5.59|) which corresponds to five European call options. The exercise price we 
take equal to ii^ = 0.914 in both cases. 

Thereafter we calculated numerically the value of the hedge-cost of the derivative 
U8{S,t) with a payoff function which corresponds to the eight European call op- 
tions with the same value E = 0.914 as before and the same expiry date T = 0.9 
and the same value p = 0.03. In these cases we use the grid with Ns = 38, 
Nf = 18, i.e., with 39 nodes in the space direction and with 18 time layers. In 
the linear case it makes no sense to calculate once more the value for this deriva- 
tive, we may just add the values U3{S,t) and U5{S,t) obtained in the former 
cases. However, in a nonlinear model where a sum of solutions is not neces- 
sarily a solution too, the difference between these two cases may be significant. 
Both, the function us{S,0) which is a solution of equation p.2|) and the sum 
^3(5, 0) + u^i^S, 0) which is not equal to any solution of equation ()1.2|) are shown 
in Figure El We expect that if in a linear case we can use linearity to compose 
solutions, in the nonlinear case we shall calculate the hedge-cost for each deriva- 
tive for its own. Indeed, in Figure El we see a strong difference between the values 
of the hedge-costs for the derivatives calculated in the linear and nonlinear cases 
in the neighbourhood of the exercise price E. 

We use the completely implicit method also for the numerical calculation of an 
hedge-cost for an option with an essential different payoff as in the case of a 
European call or a strangle. As an example we take a bull-price-spread option 
with the payoff 



u 



spreadiS, T) = max{S - El, 0) - max(S' - E^, 0), Ei < E^. (5.60) 



We used the same system of difference equations ()5.56|) - ()5.58|) and studied the 
option hedge-cost for various values of p. The results are represented in Figure 
ITUl and show the strong difference between linear and nonlinear cases of Black- 
Scholes equations as well as a strong dependence of the option hedge-cost on the 
feedback-effect for a large trader. 

All calculations were done using the program Mathematica 5.0. In order to 
solve the system of algebraic equations ()5.56|) - ()5.58|) we used the function Find- 
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Root. If we use the boundary conditions ()5.59|) . then on the interval S G [0, E] 
Ubound{S, T) = uk{S, T) = holds. Would we take as the first approximate values 
for the procedure FindRoot zeros for the values of Ui^Nt-i then this procedure 
will lead to the trivial solution for the system of equations ()5.5(i|l - ()5.58j] . To avoid 
this problem we take as the first approximate values for the procedure FindRoot 
some small constant k. We proved that the solutions to the system ()5.56|) - ()5.58p 
do not depend on this constant. In our calculations we used k = 0.03 for the 
calculations represented on Figures [71 |H1 Eland k = 1.0 for numerical solutions 
given on the Figure E3 



6 Conclusion 

We studied the symmetry properties of the nonlinear partial differential equa- 
tion p.2|l . We found the corresponding four dimensional Lie algebra ()2.1(ij) and 
the explicit representation of the Lie group ()2.20|1 - ()2.2H|1 for this equation. The 
existence of a nontrivial Lie group allowed us to obtain the invariants ()2.24p - 
fl2.25|) which can be used as new independent and dependent variables. Using 
new scaling variables we reduced the partial differential equation to the ordi- 
nary differential equation ()3.28j) . This equation possesses a solvable Lie algebra 
spanned by the infinitesimal generators ()H.29j) . Consequently we were able to 
reduce this equation to the first order differential equation fl4.34|) . We proved 
uniqueness conditions for this ordinary differential equation. We used the alge- 
bra ()3.29|) to obtain families of invariant solutions ()4.42|) and ()4.46p and proved 
that the uniqueness conditions for these solutions can fail just in the point S = 0. 
The invariant solutions have boundary conditions which approximates payoff of 
strangles. We studied sensitivity parameters for these solutions and gave graphi- 
cally representations for the dependences of these parameters on time and value 
of underlying asset. We used the obtained invariant solutions to test numerical 
methods. We proved that the best result can be obtained with a completely 
implicit method. We used this numerical method to find numerical solutions for 
calls and bull-price-spread options. In all studied cases we have seen a strong 
dependence of the option hedge-cost on the feedback-effect of the large trader. 
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Figure 1: Plot of the solution u{S,t) (IOH|l with S G (0,2], t e [0,1] and 
parameters a = 0.35, m = 0.5, p = 0.1, di = d2 = 0. 
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Figure 2: Plot of the explicit solution m(5', t), ()4.4(i|l . (dashed line) with param- 
eters p = 0.05, m = 1338.0, di = 140.0, cis = 295139, t = compared with 
the solution u{S, 0) (solid line) of the linear Black-Scholes model ()5.53p with 
the payoff Ustrangie{S,T) (thin solid line). The parameters for the strangle are 
r = 0.02, a = 0.25, Ep = 15.0, Ec = 20.0 and T = 1.0. 



The Greeks for solutions ()4.4(ij) . 
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Figure 3: Plot of pA = p-Qg with Figure 4: Plot of pT = p ^^a' with 
S G (0,100], t G [0,1] and parameters S G (0,100], t G [0,1] and parameters 
a = 0.28, m = 8.5, di = d2 = 0. a = 0.35, m = 4.9, di = d2 = 0. 
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The sensitivity parameters G and Vega for solutions ()4.46|) . 
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Figure 5: Plot of p6 = —p g^' ' with Figure 6: Plot of pVega = —p g^' ' 
S G (0, 100], t G [0, 1] and parameters with S G (0, 100], t G [0, 1] and param- 
/T — no m — —17/7, — /7„ — n ofoT-c ^ — n "i^ m — n^ri, — ri„ — n 



a = 0.2, m = -1.7, di = d2 = 0. 



eters a = 0.35, m = 0.5, di 
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Figure 7: Plot of the numerical solu- 
tion for the hedge-cost Mi(5', 0), p.2|) . 
with the payoff Ui{S,T) ()5.52|) for var- 
ious values of p. Compare the solu- 
tion ui{S, 0) with p = 0.3 (short dashed 
line), with p = 0.2 (solid line) and with 
p = 0.1 (long dashed line). 
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Figure 8: The part of the same curves 
as in Figure [7| in the neighborhood S ~ 
E. The parameters of the European 
call are a = 0.35, T = 0.9, E = 0.914, 
S = [0.1, 2]. The parameters of the grid 
are h = 0.05, r = 0.05, Ns = 38, 
Nt = 18. 
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Figure 9: Plot of the numerical solution ^8(5', 0), p.2j) . (dots) with the payoff 
Us{S,T), ()5.59p . associated with 8 European calls compared with the solution 
u{S, 0) (solid line) of the linear Black-Scholes model ()5.53|) with the same 
payoff, with the sum of numerical solutions ^3(5,0) + ^5(5', 0) (dashed line) 
of equation ()1.2|) . The parameters are S G [0.1,2.0], t G [0,T], T = 0.9, 
r = 0.02, a = 0.35, E = 0.914, p = 0.03. The parameters of the grid are 
Ns = 38, Nt = 18, T = 0.05, h = 0.05. 
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Figure 10: Plot of the hedge-cost for a bull-price-spread option with the payoff 
function Uspread{S,T), ()5.60|) . (thin solid line) for various values of p. Compare 
the solution u{S, 0) of the linear Black-Scholes model ()5.53|) (thick solid line) 
which corresponds to p = with the numerical solutions Uspread{S, 0) to the 
nonlinear equation p.2j) with p = 0.2 (short dashed line), with p = 0.1 (dots) 
and with p = 0.05 (long dashed line). The parameters for the bull-price- 
spread option are S G (20,140], t G [0,T], T = 1.0, r = 0.02, a = 0.35, the 
exercise price for the long European call is Ei = 60.0, the exercise price for the 
short European call is Eg = 80.0. For all numerical solutions the same payoff, 
volatility, expiry date and exercise prices as in linear case are chosen and the 
parameters of the grid are h = 2, t = 0.05, Ns = 60, Nt = 20. 
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